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We carefully examine common measures of dynamical heterogeneity for a model polymer melt 
and test how these scales compare with those hypothesized by the Adam and Gibbs (AG) and 
random first-order transition (RFOT) theories of relaxation in glass-forming liquids. To this end, 
we first analyze clusters of highly mobile particles, the string-like collective motion of these mobile 
particles, and clusters of relative low mobility. We show that the time scale of the high-mobility 
clusters and strings is associated with a diffusive time scale, while the low-mobility particles' time 
scale relates to a structural relaxation time. The difference of the characteristic times for the high- 
and low-mobility particles naturally explains the well-known decoupling of diffusion and structural 
relaxation time scales. Despite the inherent difference of dynamics between high- and low-mobility 
particles, we find a high degree of similarity in the geometrical structure of these particle clusters. In 
particular, we show that the fractal dimensions of these clusters are consistent with those of swollen 
branched polymers or branched polymers with screened excluded-volume interactions, corresponding 
to lattice animals and percolation clusters, respectively. In contrast, the fractal dimension of the 
strings crosses over from that of self-avoiding walks for small strings, to simple random walks for 
longer, more strongly interacting, strings, corresponding to flexible polymers with screened excluded- 
volume interactions. We examine the appropriateness of identifying the size scales of either mobile 
particle clusters or strings with the size of cooperatively rearranging regions (CRR) in the AG and 
RFOT theories. We find that the string size appears to be the most consistent measure of CRR 
for both the AG and RFOT models. Identifying strings or clusters with the "mosaic" length of the 
RFOT model relaxes the conventional assumption that the "entropic droplets" are compact. We 
also confirm the validity of the entropy formulation of the AG theory, constraining the exponent 
values of the RFOT theory. This constraint, together with the analysis of size scales, enables us to 
estimate the characteristic exponents of RFOT. 



I. Introduction 

One of the central mysteries of glass formation is the 
origin of the dramatic increase of relaxation times ap- 
proaching the glass transition temperature, Tg, which is 
commonly interpreted as an increase of the effective ac- 
tivation energy [1, 2]. Since the low temperature (T) 
activation energy typically exceeds the energy of a chem- 
ical bond, it is natural to associate this activation process 
with the reorganization of multiple atoms or molecules. 
Indeed, there is general agreement that glass-forming liq- 
uids are dynamically heterogeneous, exhibiting a signif- 
icant fraction of particles with extreme high or low mo- 
bility relative to the mean, whose positions are spatially 
correlated [3-5]. 

Even before the phenomenology of dynamical hetero- 
geneity was clearly established, Adam and Gibbs [6] (AG) 
suggested a molecular picture of this kind in 1965, along 
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with specific predictions for the relation of the configura- 
tional entropy ^conf to the relaxation dynamics. In par- 
ticular, they proposed that reorganization in a liquid oc- 
curs via hypothetical "cooperatively rearranging regions" 
(CRR), where the activation energy for relaxation is ex- 
tensive in the number of atoms or molecules that make up 
the CRR. The AG model attributes the rapid growth of 
relaxation time approaching Tg to the progressive growth 
of the CRR size on cooling. However, the AG theory 
does not include a microscopic description of the CRRs, 
or a concrete prescription for identifying them. AG fur- 
ther argued that the configurational entropy per CRR 
is roughly independent of temperature so that the CRR 
mass is inversely proportional to the configurational en- 
tropy of the fluid - a quantity that can be estimated ex- 
perimentally by the difference of the total and vibrational 
entropies. Consequently, the entropy formulation of the 
AG theory postulates that the temperature-dependent 
activation energy for relaxation is inversely proportional 
to S'conf (the "Adam-Gibbs relationship"). This model 
has proven to be highly successful to describe the T de- 
pendence of relaxation in both experiments [7-9] (where 
•S'conf is estimated from specific heat measurements) and 
computational studies [10-17] (where S'conf can be for- 
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mally evaluated from an energy landscape approach). 

The random first-order transition (RFOT) theory [18- 
20] is based upon similar ideas to rationalize the rapid 
growth of relaxation time on cooling. In particular, the 
RFOT theory formulates the problem of the relaxation 
in glass-forming liquids in terms of an "entropic droplet 
model" , or "mosaic" picture, in which the liquid is di- 
vided into metastable regions with a characteristic size ^. 
The balance between the surface and bulk free energies 
of these regions predicts a scaling relation between ^ and 
•Sconf- The overall barrier for relaxation is also assumed 
to scale with ^, providing a generalized relationship be- 
tween S'conf and relaxation. Notably, the AG relationship 
can be recovered by an appropriate limit of RFOT, so 
that these models are potentially directly linked. As in 
the case of the AG theory, RFOT theory does not provide 
a specific molecular definition of the length scale of collec- 
tive motion. Thus, both the AG and RFOT approaches 
leave the precise nature of cooperative rearrangements 
and their relation to dynamical heterogeneity open to in- 
terpretation and quantification. 

Computer simulations have been particularly helpful 
to quantify the nature of dynamical heterogeneity, as the 
spatial and temporal heterogeneity of glass- forming fluids 
is difficult to probe directly with experiments. It is now 
appreciated that atoms or molecules of extreme mobility 
(or immobility) tend to cluster, and that the most mobile 
clusters can be further divided into groups of atoms or 
molecules that move cooperatively in a roughly co-linear, 
or string-like fashion [21-27], and this phenomenon has 
been confirmed experimentally in colloidal particle track- 
ing measurements. [28-31]. Consistent with the ideas of 
the AG and RFOT theories, the sizes of clusters and 
strings grow on cooling toward Tg, but it is not clear 
if either of these structures are appropriate measures of 
the size scales envisioned by these theories. Earlier works 
have considered both the possibilities of using the mobile 
particle clusters or the strings as the CRR of AG, and 
each of these studies indicated promising results [32-35] . 
We should point out that there are other ways to char- 
acterize the length scales of heterogeneity. In particular, 
the use of a four-point correlation function offers an ap- 
proach rooted in the framework of statistical mechanics 
that reveals a growing dynamical length scale on cool- 
ing [36-44]. 

In this work, we systematically dissect cluster and 
string-like nature of the heterogenous motion in a model 
glass-forming polymer melt, and then consider what mea- 
sure or measures of dynamical heterogeneity, if any, may 
appropriately quantify the size scales envisioned by the 
AG or RFOT approaches. In doing so, we expand on a 
general methodology to identify subsets of extreme im- 
mobility. Our results span a broad temperature range, 
from very high T, to somewhat below the crossover tem- 
perature Tc often associated with mode-coupling theory. 
We find that, at the characteristic time of maximal clus- 
tering, the structures of mobile and immobile clusters 
exhibit statistical properties that are consistent with the 



properties of equilibrium branched polymers, which are 
the same as clusters approaching a percolation transition 
(lattice animals). When mobile clusters are decomposed 
into strings, the geometry of short strings are consistent 
with self-avoiding walks, while larger strings (that appear 
at low temperature) behave like simple random walks. 
We find that none of the cluster types that we study form 
compact objects when examined at their respective char- 
acteristic times. Moreover, the characteristic times of 
mobile and immobile clusters provide a physically trans- 
parent way to understand decoupling phenomena, as the 
mobile cluster time scales have essentially the same tem- 
perature dependence as diffusive time scales, while the 
immobile cluster time scale follows the structural relax- 
ation time. We consider both mobile clusters and strings 
as possible descriptions of CRR in the AG and RFOT 
models, and find that the strings - which necessarily in- 
corporate large mobility and cooperativity of displace- 
ment - best accord with the quantitative description of 
the mass or length scales of cooperative clusters described 
by both these theories. 



II. Model and Simulation Details 

Our results are primarily based on molecular dy- 
namics simulations of a melt containing 400 chains of 
"bead-spring" polymers, each chain consisting of 20 
monomers [45]. At this length, the chains are unentan- 
gied. All monomers interact via a force-shifted Lennard 
Jones (LJ) potential, truncated at 2.5ct so that dispersive 
attractions are included (cr is the LJ length parameter). 
Neighboring monomers along a chain also interact via a 
FENE spring potential to create covalent bonds. The 
FENE parameters are k = 30e and i?o = 1-5(t, chosen 
to create a mismatch in the length scale of bonded and 
non-bonded interactions, thus frustrating crystallization 
and making the model a good glass former [46]. All val- 
ues are reported in reduced LJ units. Standard units 
for temperature are recovered by multiplying T by e/fc^ , 
where ks is Boltzmann's constant. Time is given in 
units of (mcr^/e)^/^. The simulations cover the range 
of 0.3 < T < 2.5 at constant density p = 1.0. For 
all T < 1, we perform five independent simulations to 
improve statistics. Each simulation consists of an equi- 
libration run followed by data collection; the duration 
of each run is determined from the a relaxation time 
(discussed below) to ensure we sample only equilibrium 
states. Temperature is controlled via the Nose-Hoover al- 
gorithm, which is implemented via the rRESPA method 
using a time step of 0.002 for bond forces with 3 updates 
for each non-bonded force update [47]. 

The dynamics of this model (or the closely related 
model that excludes LJ attractions) have been exten- 
sively studied in previous simulations [46]. To provide 
basic characterization for subsequent detailed analysis, 
we first consider the relaxation of the coherent density- 
density correlation function F{q, t) (Fig. 1(a)). We evalu- 
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FIG. 1: Characterization of basic dynamical properties of the 
polymer melt, (a) The coherent density-density correlation 
function F{qo,t) for all T. The a-relaxation time is defined 
by F{qo,Ta) — 0.2. Symbols are shown for the lowest T to 
indicate typical intervals at which data are collected, (b) The 
non-Gaussian parameter a2 (t) at each T shows a peak due to 
the correlated motion that occurs roughly on the time scale 
of t* , defined by the maximum of 02 (i). The inset compares 
the behavior of Tq and t* . 



ate F{q, t) at the wave vector qq corresponding to the first 
peak of the structure factor where relaxation is slowest 
(except for the limit q — >■ 0). We define the a-relaxation 
time by F{qo, Tq,) = 0.2. The T-dependence of (Fig. 1 
inset) is characterized by simple Arrhenius behavior for 
T > Ta', for Tg < T < Ta, Tq grows significantly faster 
on cooling, and is well-approximated by the ubiquitous 
Vogel-Fulcher-Tamman (or Williams-Landau-Ferry) ex- 
pression 



Toexp 



DTa 



T-Tn 



(1) 



where Tq is an extrapolated divergence temperature that 
is typically slightly below the laboratory glass transition 
temperature Tg, and D characterizes the curvature (or 
fragility) of t^. From our simulations, the crossover from 
Arrhenius behavior Ta « 0.8 and Tq = 0.20 ± 0.01. For 
this system and density, the characteristic temperature 
Tc associated with power-law behavior r ~ (T — Tc)~^ 
has been estimated to be T^ = 0.35 [48]; thus we probe T 
significantly below T^. Additionally, we also know that 



[34]. For a simple polymer (like 
polystyrene) with Tg « 100 °C, the reduced units can be 
mapped to physical units relevant to real polymer materi- 
als, where the size of a chain segments a is typically about 
1 nm to 2 nm, time is measured in ps, and e « 1 kJ/mol. 

Since we will examine in detail the spatial heterogene- 
ity of the segmental dynamics, we also evaluate the non- 
Gaussian parameter a2 as a basic indicator of the time 
scale and strength of correlated motion (Fig. 1(b)). The 
peak of a2 defines the time t* , which provides a charac- 
teristic time scale of the spatially heterogeneous motion. 
The amplitude of the peak of 02 also increases, a con- 
sequence of the increasing degree of spatial correlations 
of the motion. Although it is not explicitly documented, 
it is implicit from many previous works [5, 49-52] that 
t* grows less rapidly than on cooling, as confirmed 
in the inset of fig. 1(b). In other words, these charac- 
teristic times "decouple" . This can be expected since 
t* is a diffusive time scale (see appendix Appendix A), 
and the diffusion coefficient D has long been known to 
decouple from structural relaxation. [3] AG never envi- 
sioned that glass-forming liquids should be characterized 
by multiple relaxation times; consequently, they did not 
distinguish between mass diffusion and momentum diffu- 
sion (i.e., viscous relaxation), but their language clearly 
relates to modeling mass diffusion. Fortunately, since 
these time scales maintain a fractional power-law rela- 
tion over a large range extending from Tg to Ta , the AG 
(or RFOT) approaches can be equally applied to either 
mass or momentum diffusion, a point that we expand 
upon below. 



III. Dynamical Clusters 
Glass Transition 



Approaching the 



It is widely appreciated that, below the onset tem- 
perature Ta, the dynamics become increasingly spatially 
heterogeneous approaching Tg. Regions with either en- 
hanced mobility or diminished mobility form in a spa- 
tially correlated manner, and the motions within mobile 
regions can be further dissected into more elementary 
groups that move in a string-like, cooperative fashion. In 
this section, we examine several ways to characterize cor- 
relations in mobility and analyze the geometry of these 
structures. 



A. Mobile and Immobile Clusters 

Since the distribution of particle mobilities varies con- 
tinuously, the first challenge is how to distinguish mobil- 
ity subsets. For a variety of systems [5, 51-53], it has 
been shown that choosing the subset of particles that 
that have moved farther than is expected from the Gaus- 
sian approximation at the characteristic time t* offers a 
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useful metric to identify the highly-mobile particles. De- 
pending on the system, these mobile particles typically 
account for 5 % to 7 % of the particles below T4. This 
is also true for the present system, confirmed in our own 
calculations and in ref. [52]. Accordingly, we follow the 
choice of ref. [52] where the same model was examined, 
and select mobile particles as the 6.5 % of particles with 
the greatest displacement over any chosen interval t. This 
allows us so see the evolution of mobile particle properties 
over all t, in addition to the characteristic time t* . 

At the other mobility extreme, we identify particles of 
extreme immobility. While a variety of methods to iden- 
tify low mobility particles have been explored in past lit- 
erature [5, 54-61], there is considerable variation in the 
details of these approaches. Moreover, many studies of 
immobile particles do not condition the selection on mo- 
bility, but rather on local packing considerations {e.g., 
icosahedral packing, Frank-Kasper clusters, a sufficient 
number of neighbors etc.). Such attempts are potentially 
valuable for relating structure to dynamical behavior, but 
in the present work we wish focus purely on dynamical 
considerations that should be applicable to characterizing 
mobility all glass-forming liquids, rather than any partic- 
ular fluid having its own unique type of local ordering. 
Consequently, we have developed a criterion for immobile 
particles based on the tendency for "caged particles" to 
cluster. We provide a detailed description of the method 
in appendix Appendix B to avoid breaking the flow of our 
main results. Broadly speaking, we can identify caged 
particles for any time t by those particles with displace- 
ments smaller than the (weakly T-dependent) plateau 
value observed in the mean-squared displacement. We 
then find the time at which these caged particles form the 
largest clusters and evaluate what fraction of the system 
these caged particles constitute. Similar to the approach 
for mobile-particle clusters, we fix this fraction for all t 
to track the evolution of the clustering of the immobile 
particles. Note that the fraction of immobile particles 
from this method is T-dependent, increasing from ss 5 % 
at the lowest T studied up to 11 % at T^. Above Ta, the 
cage size is not well defined. 

Having identified the most and least mobile particles 
at each interval t, we examine the average cluster size 
of the mobile {nM{t)) and immobile {ni{t)) subsets. We 
plot the cluster sizes (fig. 2) relative to the cluster size 
of the same fraction of particles chosen randomly; for 
immobile particle clusters, this eliminates the trivial T- 
dependence of immobile particle cluster size that arises 
from T-dependence of the fraction of immobile particles 
(see the appendix Appendix B for further discussion). 
We define a cluster by the group of nearest-neighbor 
particles that have a separation less than the nearest- 
neighbor distance, given by the distance of the first min- 
imum Tmin — 1.46 in the density- density pair correlation 
function. Fig. 2 shows the typical behavior for mobile 
particle clusters; namely, (riMit)) peaks at a characteris- 
tic time Im that increases on cooling, and that the peak 
value {nM{tM)) also grows on cooling, indicating an in- 
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FIG. 2: The dynamical average cluster size for (a) high- 
mobility particles {uMit)) and (b) low- mobility particles 
{ni{t)) at all T studied. The data are normalized by the 
value at t = 0, equivalent the cluster size of the same fraction 
of particles chosen at random. The definitions of the mobility 
groups are discussed in the main text. The inset of part (a) 
shows the percolation probability p of mobile particle clus- 
ters as a function of T; the dotted vertical line indicates the 
temperature where p « 0.5, a standard identifier of the per- 
colation transition in finite-sized systems [62]. For immobile 
particle clusters, p < 0.2 for all T, so that percolation is not 
prevalent. 



crease in the spatial extent of correlations. The immo- 
bile particle clusters exhibit the same qualitative trend. 
Note that we plot the average cluster size, not the weight- 
averaged cluster size; the qualitative behavior of both are 
the same. 

There are significant differences between (nAf(i)) and 
(n/(t)) to consider. First, we see that the relative peak 
size of the mobile particle clusters is larger and increases 
more rapidly on cooling than that of the least mobile 
clusters. At the lowest T studied, the mobile particle 
clusters become so large that the percolation probability 
p, defined as the fraction of configurations with a span- 
ning cluster, approaches unity (inset fig. 2(a)). If we 
define the percolation threshold by Pc — 0.5 (as is com- 
mon in finite systems [62]), the percolation temperature 
Tp w 0.32 for mobile particle clusters. Consequently, we 
likely underestimate the size of mobile particle clusters 
at the three lowest T studied. For the immobile parti- 
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FIG. 3: (a) Characteristic time where the cluster size peaks 
for mobile (black circles) and immobile (red squares) clusters. 
Note the similarity to the behaviors of Tc and t* , shown in 
fig. 1(b). Panels (b) and (c) are parametric plots that show 
that Im ~ t* and ti ~ Tc,\ the dashed lines indicate an equal- 
ity between these quantities. 



cles clusters, p < 0.2 for all T, so that finite-size effects 
should not be of concern. 

The characteristic time scales of these cluster types dif- 
fer significantly at low T. Specifically, the time scale for 
the peak of the mobile particle clusters i m is significantly 
smaller than the peak time tj of the immobile particle 
clusters (fig. 3), similar to the difference in the time scale 
between t* and t^. Indeed, parametrically plotting these 
quantities shows that Im ~ t* and tj ^ (fig. 3(b) 
and (c)), and in fact the respective quantities are nearly 
equal. The similarity between tM and t* has been previ- 
ously noted [33]. The linear scalings can be understood 
in the context of decoupling phenomena. Specifically, it 
is widely observed that Tq, grows more rapidly on cooling 
than the time scale associated with the diffusion coeffi- 
cient D, giving rise to a breakdown in the Stokes- Einstein 
relation D/T ~ - the same decoupling phenomenon 
previously discussed for t* . Decoupling can be qualita- 
tively understood as a consequence of dynamical hetero- 
geneity, since D will be dominated by the most mobile 
particles, while relaxation functions (and hence Tq.) will 
arise from the least mobile particles. Consequently, we 
anticipate, and indeed observe, decoupling between the 
ti and Im timescale that matches the decoupling between 
the t* and Tq, timescales. 



B. Cluster Size Distribution and Fractal Dimension 

We next provide a more complete account of the geo- 
metrical properties of these clusters. Figure 4 provides 
a visualization of typical mobile and immobile particle 
clusters. To quantify the structure of these clusters, we 
first consider the size distribution P{n) of the clusters 



at the characteristic times Im and t/ of the mobile and 
immobile clusters, respectively. This distribution for mo- 
bile particle clusters has been previously examined for a 
variety of systems, where it is appreciated that Pin) can 
be described by a power law with an exponential cut-off, 
namely 

P{n) ^ n'^^ exp(n/no), (2) 

where tiq is proportional to (n) . This distribution arises 
in the description of equilibrium branched polymers and 
clusters approaching a percolation transition (commonly 
referred to as "lattice animals"). We shall return to this 
analogy to help understand our findings. The Fisher ex- 
ponent Tp (using standard notation from percolation the- 
ory [62]) should not be confused with a time scale. 

Fig. 5 shows that both the mobile and immobile par- 
ticle clusters follow eq. 2, albeit with different exponents 
Tp- For mobile-particle clusters, we find Tp = 1.85 ± 0.1. 
Note that earlier work [52] for this same model indicated 
Tp 1.6, but that work was limited to much smaller clus- 
ters, and as a consequence was dominated by the behav- 
ior at small n^. Our Tp estimate is consistent with that 
for mobile particle clusters in the Kob-Andersen binary 
LJ liquid (r « 1.86) [5], and slightly larger than that for 
the Kob-Andersen lattice gas model (t sa 1.6). All these 
Tp estimates are smaller than found in percolation the- 
ory near the percolation transition in 3D {tp = 2.18) [62]. 
These variations suggest that Tp may be material depen- 
dent. The mass distribution P{n) of the least-mobile 
particle clusters exhibits similar scaling features to P{n) 
for the most-mobile clusters, but the exponent Tp dif- 
fers. In particular, Tp « 2.2 is close to that expected for 
percolation, although assignment of a precise numerical 
value to Tp for the least-mobile clusters is difficult, given 
the present data. 

We can better understand the value of Tp by consid- 
ering the possibility that mobile and immobile particle 
clusters are analogous to equilibrium branched polymers, 
which are directly related to percolation clusters. In 
three dimensions, it is known that Tp ranges from 1.5 
to about 2.2 for lattice animals [63-65] and percolation 
clusters [62, 66], respectively, so that Tp for branched 
polymers can be expected to be somewhat variable. This 
exponent reflects the effect of strong excluded volume 
interactions between and within these different types of 
model branched polymers. In mean field theory, which 
aplies above 8 and 6 dimensions for both lattice animals 
and percolation clusters, respectively, Tp is exactly equal 
its classical Flory-Stockmayer value of 5/2. Basically, 
lattice animals are swollen branched polymers and per- 
colation clusters are branched polymers with screened 
excluded volume interactions so that these structures are 
branched polymer analogs of self-avoiding and random 
walk (more precisely, 0-polymers) polymers describing 
equilibrium linear polymer chains. In short, our expo- 
nent estimates for Tp are consistent with the expected 
exponent range for branched polymers. 

To further characterize the geometrical structure of 



5 



6 Starr, Douglas, and Sastry 



Accepted for publication in Journal of Chemical Physics 




FIG. 4: Typical examples of (a) the most mobile and (b) least mobile clusters. Different clusters are shown in different colors, 
and the segments of all chains are shown translucent. Panel (c) shows the same mobile clusters (all colored red) and immobile 
clusters (all colored blue) to facilitate comparing their relative spatial distribution. 



these clusters at their characteristic peak times, we ex- 
amine the fractal dimension d/ of the clusters defined by 



the scaling of cluster size 



(3) 



where 



1 




10' 



10 10 

cluster size 



10'' 



FIG. 5: The distribution of particle cluster sizes P{n) for (a) 
mobile and (b) immobile particle clusters. The distribution 
can be described by a power-law with an exponential cut-off, 
like clusters nearing the percolation transition. The dashed 
lines indicate a characteristic power-law. Different colors rep- 
resent different T, as in previous figures. 



2N ^ 



(4) 



is the radius of gyration, and i and j denote particles 
indices within a given clusters. Earlier work has sug- 
gested that sa 2 for mobile particle clusters [67, 68], 
which corresponds to the value for lattice animals in 3D - 
that is, percolation clusters below the percolation thresh- 
old Pc- For mobile-particle clusters, we indeed find that 
smaller clusters have df ~2 (fig. 6). However, for larger 
clusters, which only occur for lower T, it appears the 
scaling crosses over to a larger 2.5. Since the ap- 

pearance of these large clusters occurs only at low T, 
extracting the best fit result for c?/ at each T results in 
the effective df growing from roughly 2 to near 2.5 on 
cooling (insets of fig. 6). For immobile clusters, the scal- 
ing of mass on Rg for small and large clusters does not 
noticeable change with size. As in the case of mobile 
clusters, df grows from 2 to near 2.5 on cooling. Thus, 
there is significant similarity in the geometrical structure 
of the mobile and immobile particles clusters. However, 
the precise values of d/ should be taken with caution, 
since the range of the data covers less than a complete 
decade in Rg. 

We can understand the changing value oi df by again 
considering the analogy to equilibrium branched poly- 
mers and lattice animals. Specifically, lattice animals 
are self-avoiding branched polymers with strong excluded 
volume interactions, and have a fractal dimension df — 2 
in three dimensions [62, 69]. The mobile and immobile 
particle clusters conform to this scaling at relatively high 
T, where they are sparse and not strongly interacting 
with each other. Like the mass distribution exponent 
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FIG. 6: Scaling of cluster mass n with the radius of gyra- 
tion Rg to define the fractal dimension rf/ of (a) mobile and 
(b) immobile particle clusters. The solid lines represent the 
simulation results, and different colors indicate different T, as 
in previous figures. Small circular symbols are shown for the 
lowest T, to indicate the density of data. The bold dashed 
lines are a guide to the eye, and provide approximate bounds 
on df. This is more clearly seen in the insets, which show the 
T dependence of df. The dashed lines of the insets indicate 
the limiting behaviors discussed in the text. 



Tp, the fractal dimension c?/ of branched polymers is 
also sensitive to excluded volume interactions. Branched 
polymers with screened, excluded volume interactions be- 
have like percolation clusters, and have a fractal dimen- 
sion df ~ 2.5 in 3D [62]. Thus, if we consider our clusters 
to be analogues of branched polymers, the increase of d/ 
upon cooling can be interpreted as the result of screening 
of the interactions. This situation is natural since, as the 
clusters grow upon cooling, the concentration describing 
the onset of their mutual interaction will decrease. We 
note that this crossover in exponent values has been an- 
ticipated for the branched polymer structures associated 
with Coniglio-Klein clusters in the Ising model [70]. 

The crossover behavior we observe for the exponent df 
with temperature is not apparent in the size-scaling ex- 
ponent Tp- This might be understood from the fact that, 
in percolation, Tp is far less sensitive than d/ to changes 
in cluster structure or dimensionality [62, 71]. Given 
that the range of data used to determine Tp and d/ is 
rather modest, further study to determine whether these 




FIG. 7: The dynamical string length L{t) for all T studied. 
The inset shows the T dependence of the characteristic peak 
value, simply denoted as L. 



clusters can be exactly identified with branched equilib- 
rium polymers is merited. Additionally, examination of 
the anisotropy of clusters will be valuable to improve the 
comparison to branched polymers, since cluster shape is 
often a more discriminating metric of the cluster type 
than the size distribution or fractal dimension [72] . 



C. String-Like Cooperative Motion 

Mobile-particle clusters can be further decomposed 
into subsets of string-like groups of cooperatively mov- 
ing particles. We now consider the properties of these 
"strings" , following an analysis parallel to that just dis- 
cussed for the clusters. To identify string formed by mo- 
bile particles, we follow the procedures originally devel- 
oped in ref. [21]. Specifically, using the same mobile par- 
ticles that we use to identify clusters, we consider two 
mobile monomers i and j to be in the same string if, 
over an interval t, one monomer has replaced the other 
within a radius S. Following ref. [22], which examined the 
same polymer model, we choose S — 0.55, although the 
results are not strongly sensitive to this choice for rea- 
sonable values of 5. Since we study a polymeric system, 
it is worth noting that the string-like collective motion is 
not strongly correlated with chain connectivity [22] , so it 
should not be confused with reptation-like motion. 

For reference, we first show the average length (num- 
ber of monomers) of a string L{t) for all T studied (fig. 7. 
As expected, L{t) has a peak at a characteristic time 
which we label i/,, and the time scale and amplitude of 
this peak grow on cooling, indicating increased coopera- 
tive motion nearing Tg. Since the strings are subsets of 
the mobile particle clusters, the peak value of L is sig- 
nificantly smaller than that of the mobile particle clus- 
ters. As a consequence, the percolation probability of the 
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FIG. 8: The characteristic time of the peak string length. 
The inset shows that, like the characteristic time of mobile 
particle clusters, tL scales nearly linearly with t* . 



strings, even at the lowest T simulated, is negligible. 

The characteristic time tj^ of the strings (fig. 8) is sim- 
ilar to Im for the mobile particle clusters, but is slightly 
larger, consistent with ref. [33]. Moreover, Hke Im, th 
scales linearly with t* (fig. 8 (b)). Since t* scales lin- 
early with the characteristic diffusion time (appendix Ap- 
pendix A), this helps to clarify that the mobile particle 
time scales captured by the clusters and strings relate to 
a diffusive relaxation time, rather than the a-relaxation 
time. This time scale is naturally shorter than as a 
consequence of the breakdown of the Stokes-Einstein re- 
lation. 

To complete the characterization of the strings, we ex- 
amine the distribution of string lengths P{L) and their 
fractal dimension at the characteristic time t^. As ex- 




String Length L 



FIG. 9: The distribution of string lengths P{L) follows an 
exponential law that can be anticipated by analogy with equi- 
librium polymerization [73]. 




Radius of Gyration 



FIG. 10: The scaling L ~ Rg^ to determine the fractal di- 
mension d/ of strings. The solid lines represent the simulation 
results, and different colors indicate different T, as in previous 
figures. Small circular symbols are shown for the lowest T, to 
indicate the density of data. The bold dashed lines are a guide 
to the eye, and provide approximate bounds on d/. This is 
more clearly seen in the insets, which show the T dependence 
oi df. The dashed lines of the inset indicate expected limiting 
behaviors. Specifically, there appears to be a modest change 
in df from « 5/3 for short strings (the approximate value for 
a self-avoiding walk in 3D) to « 2 (a simple random walk, 
also characteristic of many branched polymers). Accordingly, 
a typical string at low T axe less extended, just as in polymer 
chains in dilute versus concentrated solutions. 



pected from earlier works [21, 22], fig. 9 shows that P{L) 
follows an exponential distribution that is characteristic 
of linear equilibrium polymers [73] . To estimate the frac- 
tal dimension df, we examine the scaling between L and 
Rg in fig. 10. For short strings, we find that df 5/3, 
consistent with a self- avoiding walk in 3D. For longer 
strings, the scaling relation approaches d/ = 2, the frac- 
tal dimension of simple random walks or self-avoiding 
walks with screened, excluded volume interactions [74]. 
This screening effect has been seen in simulations of dy- 
namically associating linear chain polymers [75]. Since 
longer strings are prevalent only at low T, the effective 
df from fitting the entire range is T-dependent, growing 
from 5/3 to 2 on cooling (inset of fig. 10), reflecting an 
increased screening of excluded volume interactions upon 
cooling as in the branched dynamic clusters. Hence, the 
strings appear to become somewhat more compact on 
cooling. It has been argued that cooperative motions 
should become fully compact {i.e. df — > 3) at low T [76], 
but we see no indication of such a collapse for any of the 
cluster types we have examined. As for our data for 
mobile and immobile particle clusters, the precise values 
of df should be taken with caution, since the range of the 
data covers less than a complete decade in Rg. 

The time scale at which one examines cooperativity 
can be expected to be important in the consideration 
of the cluster geometry. Thus, we next consider the 



8 



9 Starr, Douglas, and Sastry 



Accepted for publication in Journal of Chemical Physics 




time f 

FIG. 11: The temporal evolution of fractal dimension df of 
mobile clusters and strings at the lowest T = 0.30. The 
mobile-particle cluster structure is relatively insensitive to t, 
but the string geometry differs dramatically between t* and 



fact that the geometrical structure of mobile clusters and 
strings depends on the time scale on which one examines 
these objects. Our previous analysis focused on df at the 
characteristic peak time of mobile particle clustering and 
string size, which is close to i*, a time that is significantly 
smaller (at low T) than Tq. Fig. 11 shows the temporal 
evolution of d/ for the lowest T studied. While df for 
mobile particle clusters is weakly dependent on time, df 
for the strings is indeed strongly dependent on the time 
scale considered. In fact, on time scales approaching the 
structural relaxation time, the strings appear compact 
(df ~ 3), which may explain contradictory claims that 
cooperative motions should form compact regions at low 
T [76]. On this long time scale, the strings are quite 
small. This result emphasizes the fact that it is critical 
to examine the cooperativity of motion on the appropri- 
ate time scale, and thus quantification of these scales is 
necessary. 



IV. Dynamical Scales and Relaxation 

A central challenge in describing glass formation is the 
origin of the rapidly increasing relaxation time approach- 
ing Tg . This is the defining characteristic of fragile glass- 
forming fluids. If one makes a natural assumption that 
relaxation is an activated process, transition state the- 
ory [77-79] indicates a general Arrhenius temperature 
dependence 

r = Toexp[Ai^/T], (5) 

that is often observed in condensed phase relaxation pro- 
cess and in the rates of chemical reactions. At high T, 
AF{T) = AFa, a constant, giving the widely known 




Temperature T 

FIG. 12: Temperature dependence of the activation free en- 
ergy AF{T), evaluated from the relaxation time according 
to eq. (6). The dashed line indicates the high T asymptotic 
value AFa = 2.2. The inset shows values of AF, string mass 
L, and mobile cluster mass um normalized by their values at 
Ta to facilitate comparison of the growth of these quantities 
on cooling. 



Arrhenius form. At lower T, this relation defines a gen- 
eralized T-dependent activation free energy 

AF(T) =Tlnr/ro, (6) 

which we show for our data in fig. 12. This provides a 
simple parametric description of the problem at hand: 
how can we understand an activation barrier the grows 
on cooling to a value that is several times larger than 
its high-T limit? Approaching Tg, this growth typically 
reaches 4 to 8 times the high-T limiting value AFa , and 
the exponential nature of activation leads to extremely 
large changes in relaxation. The key element to explain 
the increase of AF{T) is to recognize that such values 
cannot be readily reconciled on the basis of single par- 
ticle motion. Both the AG and RFOT approaches are 
built upon the notion that many particles are involved in 
relaxation, and the scale grows on cooling toward Tg. Ac- 
cordingly, the change of AF{T) constrains any attempts 
to explain the change in relaxation time of glass-forming 
liquids in terms of a growing dynamical size scale. Thus, 
as a simple starting point, we consider the relative growth 
of AF with those of the cluster and string sizes in the 
inset of fig. 12. We also make the mathematically trivial, 
but conceptually important point, that the existence of a 
fractional power law relation between r and t* (fractional 
Stokes-Einstein relation) implies that the reduced activa- 
tion energy applies to both relaxation times, and indeed 
all transport properties obeying such a power scaling re- 
lationship. This explains why the AG model for diffusive 
relaxation can be equally well applied to structural re- 
laxation. 
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A. Summary of the AG and RFOT Predictions 

The seminal work of Adam and Gibbs helped to es- 
tablish a picture of dynamics nearing Tg where mo- 
tion is dominated by "cooperatively-rearranging regions" 
(CRR) , thereby introducing the importance of a dynami- 
cal size scale. Both the AG and RFOT approaches build 
on the activation picture for dynamics. For most flu- 
ids, the high T [i.e. T > Ta) dependence of relaxation 
is given by eq. (5), where AF = A_Fa. AG associated 
this high T activation barrier with uncorrelated, single- 
particle motion. On cooling toward Tg, AG argued that 
motion becomes dominated by CRR, and that the barrier 
AT is extensive in the number z of rearranging monomers 
in a CRR, so that 

T~exp[zAT^/r]. (7) 

AG further argue that the fluid can be decomposed into 
N/z such CRR, each of which has a configurational en- 
tropy Sconf, SO that the total 

N 

"5'conf = Sconf- (8) 
Z 

Consequently, the relaxation can be directly relating 
•Sconf via 

T~exp[^/(r5conf)], (9) 

where the free energy A subsumes previous constants. 
This configurational entropy picture has proved highly 
successful in capturing the T dependence of many super- 
cooled fluids [7-17]. Unfortunately the CRR and ^conf 
are not explicitly defined by AG. Fortunately, numer- 
ous works have shown that a potential energy landscape- 
based definition of S'conf appears robust [10-17]; other 
studies have shown that z might be defined in terms of 
string or mobile cluster size [32-35] - a point we will 
examine in the context of both AG and RFOT. 

The similar RFOT description is built around a scal- 
ing description of the problem [18-20]. RFOT theory 
proposes a "mosaic" picture, in which the liquid is di- 
vided into metastable regions with a characteristic size ^ 
(the "mosaic" length) - conceptually like the CRR idea 
of AG. RFOT assumes the free barrier for reorganization 
has a general scaling with size, so that 

AF-e'^. (10) 

Accordingly, the implication is that t scales as 

r ~ eMe^/T] (11) 

The free energy of this "droplet" is a balance between 
the entropic contribution from the degeneracy of states 
TSc^f^ (for a compact droplet in dimension d) and sur- 
face free energy, which should scale as T(T)^^, where T 
is a generalized surface tension of the entropic droplet, 
and the suface scaling exponent 9 < d — 1. If the en- 
tropic droplet is not compact, such as is the case for our 



clusters and strings, we can generalize this argument sim- 
ply by replacing d with df, as in the application of the 
droplet models to critical phenomena. The ordinary sur- 
face tension of many fluids is often found to grow approx- 
imately linearly on cooling, and Cammarota et al. [80] 
argue that T of RFOT should grow at least linearly on 
cooling. Thus, assuming T(r) ~ T and balancing the 
surface and volume effects yields a scaling between ^ and 
configurational entropy, 

i-vsi'j^r'^. (12) 

Combining eqs. 11 and 12 yields the generalized AG-like 
relation between relaxation and entropy, 

T^eMAIiTS^or^tf/^'i-^ (13) 

We note that the concept of surface tension here is for- 
mally unclear, since there are no explicitly co-existing 
phases. However, subtle differences in the local packing 
of highly mobile and immobile particles certainly con- 
tributes an energy gradient near the interface of these re- 
gions. Similarly, refs. [81, 82] found that mobile-particle 
clusters found in a melting crystal can be identified with 
the nucleation of a fluid phase, so that the notion of a sur- 
face tension proper for mobile regions has a well defined 
meaning in this context. Further study in glass-forming 
liquids may illuminate the notion of the mobile particle 
clusters having a surface tension. Additionally, a more 
general scaling of T than simple proportionality to T 
would lead to a slightly different scaling relation between 
f and S'conf- 

The values of the exponents "0 and 9 are not fixed in 
the theory, but on general physical grounds the exponents 
should obey the inequahty 6* < -0 < d - 1 [20, 83]. The 
well-established entropy form of AG is recovered from 
RFOT provided that ip = {df — 9), leaving only one free 
exponent. The original presentation of RFOT by Kirk- 
patrick, Thirumalai, and Wolynes [18] presumes df — d 
(compact droplets), and argues that 9 = d/2 and ip — 9, 
satisfying the AG entropy form and exponent inequal- 
ity. More recently, the exponents have been examined 
in a number of computational and experimental analy- 
ses [42, 80, 84]; these studies are inconclusive regarding 
a universal value, but generally report values ~ 0-7 to 
1 and « 2 to 2.3. These values are troubling, since, 
based on the exponent inequality, we expect the surface 
scaling exponent 9 < 2 (for d — 3), and further that 9 
should be larger than ip. 

Both the AG and RFOT methods offer a way to relate 
the size or length scale of motion with relaxation times, 
but do not directly specify how this size scale should be 
measured. Hence, we now consider if the heterogene- 
ity scales defined by mobile clusters or strings might be 
appropriate for use within these theoretical descriptions. 
However, we note that, while many recent studies indeed 
focus on the size scales of heterogenous motion in this 
context, there are reasons to be cautious, and to con- 
sider other possible length scales [42, 85]. 
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FIG. 13: Confirmation that the entropy representation of the 
AG theory (eq. 9) is valid for the present system. The in- 
set shows the T dependence of ^conf and the extrapolating 
vanishing temperature Tk, which is consistent with To of the 
VFT fit (eq 1) to r. 



B. Testing the Adam and Gibbs Approach 

Although the relation between configurational entropy 
and relaxation (eq. 9) proposed by AG is not their start- 
ing point, this is the most commonly tested and most 
broadly supported prediction of the theory. Hence we 
first wish to test whether this relation is also valid in our 
system. Evaluation of Sconf is a rather cumbersome pro- 
cess, and so we describe the process fully in appendix Ap- 
pendix C. Fig. 13 verifies the validity of eq. 9. The inset 
of Fig. 13 shows that extrapolating a simple polynomial 
fit of Sconf to low T yields "Kauzmann" temperature 
Tk = 0.20 where S'conf — >■ 0. This is exactly the same 
value Tq obtained from independently fitting r by the 
VFT function. Hence, the vanishing of ^conf coincides 
with the independently extrapolated divergence of relax- 
ation time, a comforting consistency check. We note that 
Sconf is a relatively small contribution to the overall fluid 
entropy in comparison to the vibrational entropy of our 
polymer glass-forming liquid. This fact makes the experi- 
mental estimation of the difference between the total and 
vibrational contributions to the entropy particularly un- 
certain in polymer fluids, since there is no reliable means 
of estimating the vibrational entropy to high accuracy. 

We now continue to examine the proposal laid out by 
AG by considering the relation of r to a heterogeneity 
scale. As described above, the foundation of AG is that 
the activation energy for relaxation is extensive in the 
mass z of the CRR (eq. 7). Since the CRR are not de- 
fined by AG, previous works have considered whether 
the string mass L [33, 34] or the mobile particle cluster 
mass um [33, 52] might be appropriate measures. For 
water, it was shown that um has the desired behavior, 
but over a relatively limited range of r [32] ; for a simple 
spherically symmetric model, both L and rimobiic show 
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FIG. 14: Testing the mobile particles cluster and string sizes 
for consistency with the AG prediction that the activation 
energy is linear in the size z of CRR. The cluster size hm 
grows too rapidly on cooling to fit with the AG picture, the 
string size L appears consistent. See also figure 12 for a direct 
comparison of L and um to the relative change in AF with 
no free parameters in the comparison. 



the desired relation to z, but over an even more limited 
range [33]. More recently, motivated by description of 
AG that CRR are the most basic units of reorganization, 
refs. [34, 35] found L is an appropriate measure of z, but 
did not consider cluster size um- Unfortunately, none of 
these works could definitively exclude other measures for 
z. It has also been appreciated that other length scales 
associated with heterogeneity, such as from a four-point 
density correlation function, would be too large at Tg to 
be consistent with z [43]. We shall return to this point 
in the conclusion. 

Here, we check the plausibility of both cluster and 
string size as a measure of z over a substantially broader 
range of r than previous works to provide improved clar- 
ity. We thus consider substituting for z the peak string 
size L/L{Ta) or the peak mobile particle cluster size 
nM /nuiTA), where we normalize by the value at Ta so 
that z has the expected value near unity for T > Ta- 
Figure. 14 shows that, for T ^Ta logr is linear if we use 
i as a proxy for z, but not using um- The deviation from 
cluster size is a consequence of the fact that the mobile 
particle cluster size grows noticeably more rapidly than 
the effective activation free energy, a fact already appre- 
ciated in the inset of fig. 12. If AG is assumed correct at 
high T, such an exponential relation should continue into 
the T > Ta range, where z should be near unity, or at 
least reach its asymptotic value. In such a case, the value 
of AF from the low T fit using eq. 7 should be compa- 
rable to the high T estimate of AFa from an Arrhenius 
fit. Using z — L/L{Ta), we estimate AF — 1.8, some- 
what smaller that the value estimated from the Arrhenius 
fit, AFa = 2.2. One interpretation of this discrepancy 
is that using L{Ta) as the normalizing factor is not en- 
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FIG. 15: Testing Sconf oc 1/L. The inset checks for the best 
power-law relation, which yields an exponent slightly larger 
than one. 
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FIG. 16: The mean radius of gyration (Rg) for the strings and 
mobile particle clusters. This defines a characteristic length 
that we can test within the RFOT framework. 



tirely correct, and a somewhat smaller value would be 
more appropriate. 

Given the apparent success of the string mass L to 
describe r and the broadly reported validity of the en- 
tropy formulation of AG, we test for consistency between 
these representations by checking the expected relation 
5'conf cc 1/z (fig. 15). The data support the linearity of 
the relationship, although there are some systematic de- 
viations at both the lowest and highest T. This suggests 
that, while L captures the generally expected behavior 
of the CRR, a more detailed refinement of the determi- 
nation of L may provide a more accurate description of 
CRR. 

The string length L appears to be the most quanti- 
tatively valid choice for the CRR, and is also consistent 
with the qualitative philosophy of AG. Specifically, re- 
call that AG envision the CRR are dominated by the 
smallest group of cooperatively moving monomers. The 
physical motivation for this is that probability of such a 
group will diminish exponentially with the size, so that 
the smallest possible group that allows for rearrangement 
will dominate the relaxation. The strings are both the 
smallest such unit, and also the only candidate in which 
all particles move in a cooperative fashion. While the 
particles of mobile particles clusters are obviously spa- 
tially correlated, there is no a priori cooperativity in 
their displacements. The strings are precisely the man- 
ifestation of mobile particle cooperativity. However, we 
should be careful to point out that, while the qualitative 
language of AG is appealing, in the end the quantitative 
predictability is the most important measure. We shall 
next explore which measure best quantitatively fits with 
the formulation of the RFOT theory. 



C. Testing the RFOT Approach 

As discussed in the previous subsection, the validity of 
the entropy formulation of AG dictates that i/j = [dj — O). 
Hence there is only one free exponent in the RFOT for- 
mulation. We shall consider two approaches to determine 
these exponents, which provides an internal consistency 
check. 

A simple, but significant, difference between AG and 
RFOT is that RFOT refers to a length scale of coopera- 
tive motion, rather than an extent or mass of collective 
motion. Consequently, to test whether any of the clus- 
ters or strings might be appropriate, we need to consider 
a length scale that defines the size of mobile clusters or 
strings. The natural length scale for these objects is the 
radius of gyration Rg at their respective characteristic 
times, which we show in fig. 16. Hence we can directly 
evaluate the exponent ijj from the scaling of t with {Rg) 
(eq. 11) for the strings and clusters. 

Figure 17 shows the scaling of r with {Rg) for the 
strings and clusters, from which we obtain with the best 
fit for the exponent ij:. Given our previous findings for 
the AG approach indicating that string mass relates to 
T while the cluster mass does not, we would expect a 
superior fit for the strings. Instead, we find that the 
exponent ij] « 1.3 for strings and clusters is identical 
within the limits of our determination. Essentially, this 
is a consequence of the fact that the T-dependence of Rg 
is nearly the same for strings and for clusters, while the 
T-dependence of cluster mass differs noticeably. This ap- 
parent paradox can be resolved by recognizing that the 
largest dimension of a cluster dominates Rg , so that sys- 
tems with different mass can have similar Rg. This is 
consistent with the expectation that the largest dimen- 
sion of mobile clusters is associated with long, string-like 
cooperative regions. 

We next check for consistency of the value with the 



12 



13 Starr, Douglas, and Sastry 



Accepted for publication in Journal of Chemical Physics 




5 10 15 

A(Ry'/T 



FIG. 17: Testing eq. 11 to determine the exponent i/) of 
RFOT. Both the mobile particle clusters and strings yield 
a consistent fit with tp « 1.3. The data for the mobile par- 
ticle clusters are shifted by two units along the abscissa for 
clarity of the figure. 



expectation that ip = df — 9. We can independently de- 
termine df — 9 from the scaling of S'conf with Rg (eq. 12), 
as shown in fig. 18. While the data deviate from a power- 
law at high T, the lower T data indicate df — ~ 1.3, 
consistent with our estimates of tj}- The success of these 
independent approaches significantly increases our confi- 
dence in these estimates. Based on our previous findings 
for df for strings, we can also estimate 9 — 0.3 to 0.7. 
The value of 9 is small in comparison with values esti- 
mated refs. [42, 80, 84]. However, a small value for 9 
is physically plausible. For example, in the Ising spin 
glass, a direct evaluation yields 9 « 0.2 to 0.35 in three 
dimensions [83, 86, 87]. Moreover, our value obeys the 
expected inequality 9 < < d — 1, which the previous 
estimates violate [42, 80, 84]. 

While both string and mobile particle cluster sizes 
demonstrate reasonably scaling within RFOT, it appears 
the success of the cluster description is dependent on 
the limiting dimension dictated by the string size. The 
strings also appeared to be the only reasonable descrip- 
tion of CRR within the AG framework. Hence, we can 
find a satisfying consistency for both the AG and RFOT 
descriptions using the strings as a measure of CRR or 
mosaic scale, where the exponent values of RFOT are 
constrained to satisfy the formulation of AG. 



V. Discussion and Conclusion 

We have examined the geometrical structure of clus- 
ters and string-like cooperative motions in a model glass- 
forming polymer melt. We found an aesthetically pleas- 
ing symmetry in the geometry of high- and low-mobility 
clusters, i. e., they both conform to statistical geometry 
of equilibrium branched polymers. In doing so, we also 




FIG. 18: Evaluation of the surface scaling exponent 6 from 
the scaling law eq. (12). The value df — 9 ^ 1.3 is consistent 
with ?/', and with the expectations from AG. 



developed a novel method to identify low-mobility par- 
ticles based on persistent caging. Most importantly, we 
have examined the question of whether these heterogene- 
ity scales can be identified with the scales anticipated by 
the AG and RFOT descriptions of glass formation. We 
found the strings apparently provide the most consistent 
description of the CRR or mosaic length described by 
these respective theories. 

An important observation arising from our work is that 
these different quantifications of heterogeneous dynam- 
ics in fact correspond to distinct relaxation time scales, 
and therefore distinct processes of importance in the re- 
laxation of a glass former. In other words, there is no 
single or unique heterogeneous dynamical scale in the 
system. These immobile particles are apparently related 
with the breakdown of the linear scaling between diffu- 
sive and viscous relaxation, while the fragility of glass- 
formation is apparently more related to string-like coop- 
erative motion [34]. Thus, conventional wisdom regard- 
ing the role of heterogeneous dynamics on typical aspects 
of glassy behavior, requires further examination. Appar- 
ently, there is no single dynamic heterogeneity scale in 
glass-forming liquids. 

There is further evidence of these scales and their sig- 
nificance on condensed matter relaxation. For example, 
recent simulations of superheated Ni crystals [82] also 
find a large increase of the non-Gaussian parameter, mo- 
bile particle clusters, and string-like collective motion. 
However, in this system, there are no immobile particle 
clusters of finite extent present; there are only immobile 
atoms in a crystal lattice, and mobile particles having the 
usual constituent strings. Significantly, there is no decou- 
pling of structural relaxation from the self-intermediate 
scattering function and the diffusion coefficient, nor any 
stretched exponential decay of the self-intermediate scat- 
tering function in the superheated crystal. Similarly, in 
a recent study [85, 88] on the dependence of dynamics 
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on spatial dimensionality, it is found that the degree of 
Stokes-Einstein breakdown decreases while the fragility 
paradoxically increases with spatial dimensionality for 
dimensionality greater than two. The findings of these 
studies suggest that the immobile particle clusters, rather 
than mobile particles, are primarily responsible for de- 
coupling and stretched exponential stress relaxation of 
glass-forming liquids. This possibility merits systematic 
study and points to the different types of heterogene- 
ity (mobile and immobile particles) having significantly 
different impacts on the fluid dynamics. In other words, 
dynamic heterogeneity comes in different types that must 
be properly discriminated. 

As we alluded to earlier, another common approach to 
extract a length scale for heterogeneity is via a four-point 
correlation function. Proper determination of the four- 
point scale ^4 can be strongly affected by finite size [43], 
but careful extractions have show that ^4 grows more 
rapidly than would be expected for the CRR of the AG 
theory [42] . It is possible that ^4 could be consistent with 
the mosaic scale of the RFOT theory, but we expect a sin- 
gle measure should be compatible with both approaches, 
since they are largely complementary. The reason for the 
difference in the scaling of ^4 with that observed for the 
strings can be readily understood by considering their 
characteristic times; ^4 is determined at the time of the 
peak in the four-point susceptibility, which has essen- 
tially the same temperature dependence as r^. As we 
have shown, Tq. is also the time scale of immobile parti- 
cles, and is distinct from the time scale (and hence length 
scale) of mobile particles and strings. Accordingly, ^4 is 
primarily sensitive to particles of low mobility [59] . This 
is a consequence of the fact that, in defining the four- 
point function, a particle size a is introduced to limit the 
effects of vibrational motion, and is chosen to be larger 
than the typical cage size, following ref. [39]. The choice 
of a controls the scale of relaxation associated with the 
four-point function. Choosing a smaller value of a, closer 
to the cage size, for example, should lead to a measure 
of heterogeneity on a smaller time scale, perhaps similar 
to t*; such a choice would presumably be more sensi- 
tive to string-like excitations. Efforts in this direction, 
along with other approaches to extract the size scale of 
string-like cooperativity, would be valuable to better un- 
derstand the findings of the present paper within a more 
traditional liquid-state correlation function approach. 

In conclusion, the analysis of the relationship between 
the various clusters in the context of the AG and RFOT 
theories reveals that the strings are a particularly good 
candidate for the CRR of AG theory. While a similar 
conclusion may also be reached in the context of RFOT, 
the fact that all the heterogeneous clusters considered 
here show fractal structure with a fractal dimension lower 
than d = 3 makes a conclusive comparison difficult at this 
time, since the mosaic picture in the RFOT framework 
normally assumes the rearranging regions to be compact. 
Consequently, it will be valuable to revisit the formula- 
tion of RFOT, as the suggested fractal nature of the en- 



tropic droplets has implications for the concept of the 
effective surface tension of these regions. 



Appendix A. Non-Gaussian and Diffusive 
Time Scales 

The non-Gaussian parameter 012 (0 is often used to 
quantify the deviation of particle or segmental displace- 
ments from the Gaussian distribution expected for simple 
fluids. The maximum deviation occurs at a characteris- 
tic time t* , and it is well known that t* is smaller than 
Ta, and has a weaker temperature dependence than r^. 
It is then natural to wonder to what physical process t* 
relates. 

Combining the Stokes-Einstein relation for spheres 



D 



(Al) 



where D is the diffusion coefficient, rj is the fluid vis- 
cosity, and Rh the particle hydrodynamic radius, with 
Maxwell's relation r = ?y/Goo (where Goo is the high- 
frequency shear modulus) leads to 



D 1 

7^ c>c — . 

T T 



(A2) 



In other words, the reduced diffusion coefficient D/T 
should define an inverse relaxation time. Given the rel- 
ative slow variation of t* with T compared to from 
the intermediate scattering function, we check whether 
t* can be identified with a diffusive relaxation time de- 
fined in this way. Due to the polymeric nature of our sys- 
tem, D cannot be readily evaluated, since mean-squared 
displacement will only be linear on a much longer time 
scale, associated with the chain center-of-mass diffusion; 
this requires chain displacements at least on the order of 
the chain radius of gyration, more than we can readily 
simulate at low T . However, we can check for a rela- 
tion between D/T and t* for the Kob- Anderson binary 
Lennard- Jones fiuid [49], the most commonly studied 
computational glass-forming system, where low T data 
is accessible. 

Figure 19 shows that the (inverse) characteristic dif- 
fusion time D/T is linear with t* for the entire range of 
data, covering several decades. This result clarifies that 
t* can be associated with a diffusive time scale. Given 
the known "decoupling" of structural relaxation and D/T 
from the Stokes-Einstein relation, we accordingly expect 
the same decoupling between t* and in our polymer 
system, as observed in fig. 1. 

It should be appreciated that the decoupling relation 
between D for the overall chain displacements of a poly- 
mer and the long-time shear-stress relaxation time can 
exhibit a separate relationship from t* and Tq, [89] . This 
is a consequence of the fact that heterogeneity at the 
scale of the chain radius of gyration can differ from het- 
erogeneity at the monomer or segmental scale. Thus, in 
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O Kob-Anderson Binary LJ system 

— r ~ (D/r)-' 



10 

Dl T 




10 10 10 



10' 10 10 

time 



FIG. 19: Linear scaling between t* and the (inverse) char- 
acteristic diffusion time D/T in the Kob-Anderson binary LJ 
mixture [49] , demonstrating that the time scale t* corresponds 
to the time scale for mass diffusion. 

the polymer system, t* should be thought of as relating to 
a local monomer diffusive process, and Tq to a segmental 
structural relaxation time. 



Appendix B. Immobile Particle Definition 

To study the structure of highly immobile particles, we 
need to devise a physically sensible algoritm that picks 
out an appropriate subset of low-mobility particles. Since 
the caging of particles by their neighbors is one of the 
hallmarks of glass formation, we utilize the concept of 
"caged particles". To do so, we must identify the cage 
size. We can formally do this via the mean-squared dis- 
placement (r^(i)). Fig. 20(a) shows that (r^(t)) has a 
plateau at a characteristic size in the approximate range 
0.04 to 0.09 (for a cage radius 0.2 to 0.3). To precisely 
define the cage size, we take advantage of the fact that 
the logarithmic derivative d{\n{r'^{t)))/d{\nt) exhibits a 
clear minimum on the time scale of particle caging, tcagc- 
We thus define the cage size by Tcagc = (f^(icage))"^^^- 
We show the T dependence of rcagc in fig. 20(b) for 
T ^ Tj^; aX higher T, {r'^{t)) transitions from ballistic 
motion (d(ln(r^(t)))/d(lni) = 2) to sub-diffusive motion 
(d(ln(r^(t)))/(i(lni) « 0.6 without intervening particle 
caging. The sub-diffusive behavior is well known for this 
model [46], arising from polymeric effects. 

Having unambiguously defined a cage-size, we proceed 
to track the behavior of caged particles (z. e. particles 
with displacement less than rcage)- Fig- 21 shows the 
fraction of caged particles, which, as expected, decreases 
with time. We note that, formally, this fraction is iden- 
tical to the "self-overlap" Qs{t) used in the four-point 
correlation function formalism [39] , although the particle 
"size" used is typically fixed at a value 0.3, independent 



T 




FIG. 20; (a) Identification of the cage size from the mean- 
squared displacement. The inset shows the logarithmic 
derivative, which shows a minimum that we use to define 
('"cage)- (b) The resulting T dependence of the cage size 
i^'cugc}^^^ ■ Particles with a squared displacement less than 
('"cage) are defined as "caged particles" . 



of T, and substantially larger than the cage size at low 
T. 

Since we wish to understand the tendency for these 
caged particles to be spatially correlated, we evaluate the 
average cluster size of these immobile particles. However, 
we must take into account the fact that the number of 
these caged particles decreases with time, and thus there 
is a trivial effect on the cluster size of the number of caged 
particles. To remove this trivial effect, we normalize the 
caged particle cluster size by the cluster size of the same 
fraction of particles chosen at random; this allows us to 
see how the tendency to cluster compares to the random 
case, independent of the nature of the underlying dynam- 
ics. 

Figure 22 shows the normalized cluster size of the 
caged particles as a function of time. The qualitative 
behavior matches what we observe for other dynamical 
clusters: as smaller times, the effect of clustering is weak. 
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FIG. 21: Dynamical fraction of caged particles for T < Ta- 



and there is a characteristic time where the cluster size 
reaches a peak. As discussed in the main text, this char- 
acteristic time is similar to the a-relaxation time. This 
time, along with the peak size of the caged-clusters, de- 
fine the characteristic features that we wish to capture. 

To simplify the analysis of immobile particles and draw 
a parallel to the analysis of the mobile particles (where 
there is a fixed, f-independent fraction of particles consid- 
ered), we consider a simplification of this approach that 
still captures the characteristic peak time and amplitude 
of caged particles. Specifically, if we look at the charac- 
teristic peak time of the caged particles, we can identify 
the fraction of caged particles at this time, which we 
show in fig. 23. This characteristic fraction of caged par- 
ticles increases with T, as the mobility subsets become 
less distinct at higher T. For all subsequent analysis, 
we use this T-dependent fraction for all time, to parallel 
the approach for the mobile particles. By construction. 



2.5 r- 
2.25 - 

N 

CO 2 - 



10° 10' 10^ 10^ 10* 10^ 

time 



FIG. 22: The normalized cluster size for caged particles for 
T < Ta, where the cage is well-defined. 



FIG. 23: The characteristic fraction of low mobility particles 
for several different approaches. All approaches yield similar 
fractions for all T. 



this fixed fraction reproduces the characteristic time and 
size the our time-dependent fraction of caged-particles 
reveals, so that we do not alter these important features. 

A natural concern is the sensitivity of this approach 
to the definition of the cage size, since this is the only 
parameter that must somehow be chosen. To test this, 
we also considered a T-independent cage size rcage of 0.2 
or 0.25. We find that these T independent sizes of course 
yield quantitatively different results, but that the time 
and size scales of the immobile particles all scale in the 
same way. Additionally, as shown in fig. 23, the charac- 
teristic fraction of caged particles is nearly the same as 
our definition based on a T-dependent cage size, showing 
the precise definition of cage size does not strongly affect 
the characteristic immobile fraction. 

We also considered another approach to extract the low 
mobility subset. This approach is motivated by ref. [52], 
where they determined a characteristic fraction of highly 
mobile particles by finding the fraction that maximizes 
the cluster size formed by those highly mobile particles 
relative to same fraction of particles chosen at random. 
This measure provides a way to capture the fraction of 
particles that most strongly exhibit spatial clustering, 
without consideration of the underlying mobility distri- 
bution. We use the same approach, but finding the frac- 
tion of least mobile particles that maximizes the relative 
clustering. We show this characteristic fraction together 
with our estimates from the cage size in fig. 23, and 
find that we recover nearly the same fraction by this ap- 
proach. Apparently, the characteristic fraction of immo- 
bile particles is not strongly sensitive to the exact choice 
of parameters, providing confidence in the robustness of 
our analysis. For the calculations presented in the main 
body of the manuscript, we simply use the mean of all 
these estimates. 

Finally, we point out that in all approaches to identify 
low mobility particles, there appears to be an unantici- 



16 



17 Starr, Douglas, and Sastry 



Accepted for publication in Journal of Chemical Physics 



pated feature at T w 0.35, previously estimated as the 
characteristic for this system at this density [48] . This 
is a curious result, but at present we have no explanation 
for such behavior. 



Appendix C. Calculation of Configurational 
Entropy 

Our goal is the evaluation of the configurational en- 
tropy S'conf J which enumerates the density of stable po- 
tential energy minima sampled by the melt at equilib- 
rium. Procedures for evaluating S'conf have been devel- 
oped applied to a variety of systems, including water [15], 
binary LJ mixtures [10, 11], silica [16], and orthoter- 
phenyl [17]. We follow a similar procedure, whereby the 
overall entropy can be partitioned into vibrational Svib 
and configurational components, i.e. 



5 = 5, 



conf 



S-xr 



(CI) 



Our approach will be to evaluate directly S and Svib, and 
by their difference S'conf- 



1. Total Entropy 

To evaluate the absolute entropy of the polymer, we 
employ the thermodynamic integration technique [47] . In 
this method, the free energy is calculated by parametri- 
cally coupling the potential energy of the system U^oiy to 
the potential energy U-^ci of a reference system for which 
the free energy can be directly, analytically evaluated. 
The coupling potential is normally of the form 



[/(A) = (1 - A)C/polymcr + AC/ref 



where the coupling parameter < A < 1- The free energy 
can then be evaluated by 



poly 



rcf 



(f/rcf - C^poly) d\. 



(C3) 



This procedure is complicated by the FENE poten- 
tial that bonds nearest neighbors, since it diverges as the 
bond length approaches Rq. The normally chosen ref- 
erence potentials do not restrain the bond length, and 
so the contribution from the FENE potential diverges 
strongly as A — > 1. The avoid this complication, we per- 
form a "two-step" thermodynamic integration, where we 
first perform an integration from the FENE potential to 
a harmonic bond potential 



^harm , 



Ri) 



(C4) 



that does not exhibit the strong divergence as bond 
length grows; when then perform an integration from 
the harmonically bonded polymer to the reference sys- 
tem. We choose /charm = 980 and i?i = 0.87 so that the 



location of the minimum and the curvature at the min- 
imum are very near to that of the FENE bond (when 
combined with the core LJ repulsion). While the free 
energy for the harmonically bonded polymer is not ana- 
lytically known, it is not needed, as it drops out in the 
final expression for the free energy 



^poly — -Fi-of — I {U^ 

Jo 



polymer ^harm, 
1 



dAi 



(C5) 



(C/rct - C/l larm 







For the reference potential, we use a potential that 
shares the system periodicity [90], namely 



f/rcf(r) 



3 

-Uq ^ COS 

i=l 



27r 



(C6) 



where Uq = 10 is the amplitude of the potential, L is the 
length of the container and is the coordinate of a par- 
ticle in direction i. For an N particle system interacting 
through Ureij evaluation of the partition function shows 
that 



-Fref — 



3N 



In 



pa A 



/? lIo(/?Uo) 



(C7) 



Here p is the number density and Iq (x) is the modi- 
fied Bessel function of the first kind. Combining with 
Eq. (C5), we have F for the our system at some fixed 
temperature Tq and density. Accordingly, we can evalu- 
ate the entropy for a reference To, 



EjTo) - FjTo) 
To 



(C2) We obtain S for any T by exploiting the fact that 



(C8) 



(C9) 



so that 



S(T) 



f \ dT ), 



To 



dT + SiTo). 



(CIO) 



The integrand can be evaluated numerically from data 
for E(T). Since we must explicitly include Planck's con- 
stant ?i, we select units appropriate for a monomer of a 
typical polymer, like polystyrene; specifically, we choose 
e = 1 kJ/mol, (7 = 1 nm, and m = 100 g/mol; using 
these units, Ti = 0.0635078 kJ-ps/mol. 



2. Vibrational Entropy 

The vibrational component of the entropy Svib reflects 
the contributions of the basin shape to the vibrational 
behavior. We can partition 



Sv 



Sh arm ~l" Sanh 



(Cll) 
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FIG. 24: (a) The total entropy S and vibrational entropy as a 
function of temperature, (b) The configurational component 
Sconi = S — Svib- The line is a guide for the eye. 



In the harmonic approximation, the basin entropy 



'5'harm — 



3Af-3 



1 - In 



(C12) 



where a;„ are the normal modes associated with the basin 
minimum. To evaluate ujn, we must first determine the 
basins associated with the equilibrium liquid. To do so, 
we perform a conjugate gradient minimization of the po- 
tential energy, starting from instantaneous snapshots of 
the equilibrium polymer; we locate the corresponding 
minimum, or inherent structure (IS), within a numeri- 



cal tolerance of 10 Using the configuration at the 
minimum, we evaluate the Hessian matrix 



(C13) 



the matrix of the curvatures of the potential energy. In 
the harmonic picture, the eigenvalues A„ = mw^, so we 
directly obtain the normal modes {ujn}- For each T, we 
generate the IS and {w„} for at least 100 configurations 
that are well-separated temporally. 

The anharmonic contribution to S'vib for many systems 
is negligible. To check the anharmonic contribution 5anh, 
we first consider and anharmonic energy 



E^^^{T)=E - -NkeT-eis 



(CM) 



where eis is the inherent structure energy and 3/2 NksT 
is the contribution for a harmonic solid. We can then 
evaluate 



Sanh{T) 



1 dE, 



aiih 



T dT 



dT. 



(C15) 



To obtain a valid estimate of -Banh for a basin, we must 
heat the IS very rapidly to insure that the system can- 
not change basins while heating. We perform such heat- 
ing for at least 100 IS generated from initial equilibrium 
configurations at T = 0.31. In principle, -Eanh (and thus 
•Sanh) depend on the equilibrium T from which the IS 
are obtained; in practice, we find that the IS from differ- 
ent T have a nearly identical density of states p{uj) and 
^'anh {T) . Hence we can use the behavior of i?anh {T) from 
one set of IS for any T. We find that the contribution an- 
harmonic contribution is rather small and negative, and 
can be well described by i^anh = -0.065 T for T < 0.8; 
the results in 5a„h(r) = -0.065 InT. 
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